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ABSTRACT 

Context. It is established that both radiative transfer and magnetic field have a strong impact on the collapse and the fragmentation of 
prestellar dense cores, but no consistent calculation exists yet at such scales. 

Aims. We perform the first Radiation-Magneto-HydroDynamics numerical calculations at a prestellar core scale. 
Mefhods. We present original AMR calculations including magnetic field (in the ideal MHD limit) and radiative transfer, within the 
Flux Limited Dilfusion approximation, of the collapse of a 1 M© dense core. We compare the results with calculations performed with 
a barotropic EOS. 

Results. We show that radiative transfer has an important impact on the collapse and the fragmentation, through the cooling or heating 
of the gas, and is complementary of the magnetic field. A larger field yields a stronger magnetic braking, increasing the accretion rate 
and thus the effect of the radiative feedback. Even for a strongly magnetized core, where the dynamics of the collapse is dominated by 
the magnetic field, radiative transfer is crucial to determine the temperature and optical depth distributions, two potentially accessible 
observational diagnostics. A barotropic EOS cannot account for realistic fragmentation. The diffusivity of the numerical scheme, 
however, is found to strongly affect the output of the collapse, leading eventually to spurious fragmentation. 

Conclusions. Both radiative transfer and magnetic field must be included in numerical calculations of star formation to obtain realistic 
collapse configurations and observable signatures. Nevertheless, the numerical resolution and the robustness of the solver are of prime 
importance to obtain reliable results. When using an accurate solver, the fragmentation is found to always remain inhibited by the 
magnetic field, at least in the ideal MHD limit, even when radiative transfer is included. 
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1. Introduction 

Understanding star formation is one of the most challeng- 
ing problems in contemporary astrophysics and numerical cal- 
culations provide a useful approach to investigate it thoroughly. 
Thanks to the steadily increasing computer performances, nu- 
merical calculations can integrate more and more physical pro- 
cesses. Among these latter, the coupling between matter and ra- 
diation is a major issue. Up to the formation of the first Larson 
core (Larson 1969) l, the accreting gas can freely radiate into 
space and is nearly isothermal (optically thin regime). Once the 
gas becomes dense enough (p > lO"'-' g cm"-'), the radiation 
is trapped and the gas starts to heat up (optically thick regime). 
The transition between these two regimes controls the collapse 
and the fragmentation of the cloud. The cooler the gas, the more 
important the fragmentation. Another key issue in star forma- 
tion is the role of the magnetic field, as dense cores are observed 
to exhibit coherent magnetic structures (e.g. |Heiles & Crutcher| 
|2005] l. 

Several authors have investigated the impact of radiation 
(|Boss et al.||2000t [WhTtehouse & B ate 2006; |A ttwood et"aL] 
|2009| IKrumholz et al.||200 7; Bate 2009; Offner et al.||2009|l 



|Pu dritz|2006l|Hennebelle & Teyssier|2008||Machida et al.|2008 



and of magnetic field (e.g. [Mellon & Li, ,2009; |Banerjee & 



Pric e & Bate||2007| |2008|l on the star formation process. So 
far, however, only Price & Bate ( 2009)1 have performed calcu- 



lations with both magnetic field and radiative feedback, with 
an SPH method. Yet, these calculations explore relatively large 
scales. In this Letter, we present the first full Radiation-Magneto- 
HydroDynamics (RMHD) calculations of the protostellar col- 
lapse at small scales. We investigate in detail the impact of both 
magnetic field (in the limit of ideal MHD) and radiative feed- 
back on the fragmentation process and on the launching of out- 
flows. The results are compared with the ones obtained with a 
barotropic equation of state (EOS), that crudely mimics the tran- 
sition from the isothermal to adiabatic regime, in order to assess 
the importance of a proper treatment of radiative transfer on the 
collapse. 

2. Numerical method and initial conditions 



Send offprint requests to: B. Conmierfon 



We use the RAMSES code ( |Teyssier| [2002) 1 based on a 
Eulerian formalism (grid based method with Adaptive Mesh 
Refinement). We solve the transfer equations in the Flux Limited 
Diffusion (FLD, Minerbo 1978 ) in the comoving frame (fluid 
frame) to evaluate the radiative quantities. The ideal MHD 
equations are integrated using an unsplit second-order Godunov 
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Fig. 1. Case yu - 20: Density and temperature maps in the equatorial plane at time t= 38.6 kyr: FLD approximation with LF (left), 
FLD with HLLD (middle) and with the barotropic EOS (right). 



scheme (| Fromang et al.|2006| l. Coupling terms between matter 
and radiation, as well as radiation transport, are integrated im- 
plicitly in order to handle the very short heating, cooling and 
diffusion time scales. Details on the numerical method and its 
implementation will be presented in a forthcoming paper. 

Calculations have been performed using either the rather dif- 
fusive Lax Friedrich (LF) Riemann solver or the more accurate 
HLLD Riemann solve r (Miyoshi & Kusano]p005 ). Following 
up on former studies (Commergon et al. 20081, we impose at 
least 10 cells per Jeans length as a grid refinement criterion 
(parameter A^j). The initial resolution of the grid contains 64^ 
grids. We use the low temperature opacities of Semenov et al. 
( |2003 |l, parametrized as functions of the gas temperature and 
density. In contrast, the non radiative calculations have been per- 
formed with a barotropic EOS P/p - c^()[l + (p/Pad)'^ where 
Pad - 2.3 X 10"^^ g cm~^ is the critical density at which the gas 
becomes adiabatic, -y (set up to 5/3) is the adiabatic index and Cso 
the isothermal sound speed. 



We adopt initial conditions, similar to those chosen in pre- 
vious studies (e.g. Commer9on et al. 2008| l. We consider a 
uniform-density sphere of molecular gas, rotating about the z- 
axis with a uniform angular velocity. In the present study, the 
prestellar core mass is fixed at Mo = 1 M0 and the temperature 
at 11 K, which corresponds yo Cso ~ 0.19 km s"'. To promote 
fragmentation, we use an m = 2 azimuthal density perturbation 
with an amplitude of 10%. The magnetic field is initially uniform 
and parallel to the rotation axis. The strength of the magnetic 
field is expressed in terms of the mass-to-flux over critical mass- 
to-flux ratio fi = (Mo /<!>)/ (Mo /<i>)c. The initial energy balance is 
determined by two dimensionless parameters, namely the ratio 
of the thermal over gravitational energies a - 0.37, and rota- 
tional over gravitational energies /3 - 0.045. The corresponding 
free-fall time is % ~ 33 kyr. 



3. Results 

3.1. Casefi^ 20 

We first present the results for ju - 20, i.e. a weakly mag- 
netised core. We compare FLD calculations performed with the 
two Riemann solvers (LF or HLLD) with calculations performed 
with the barotropic EOS and the LF solver. The grid refinement 
criterion is A^j = 15. Performing similar barotropic calculations, 
Hennebelle & Teyssier ( |2008) do not report any fragmentation. 

Figure [T] portrays temperature (top) and density (bottom) 
maps in the equatorial plane for the three aforementioned calcu- 
lations at the same time f ~ 38 kyr. The black contours represent 
the transition between the optically thin and thick regions (thin 
contour) and the quasi adiabatic regions. The FLD case with the 
LF solver yields a multiple fragmentation, with a central frag- 
ment of mass 5.2 x 10"^ Mq and several (depending on the time 
and resolution) orbiting fragments of mass ~ 2 - 4 x lO"-' M© 
with separations ranging from ~ 40 AU to ~ 100 AU, whereas 
no fragmentation occurs in the two other cases. The orbiting 
fragments in the FLD-LF simulation are quite warm (~ 40 K) 
while the disk is cold (~ 1 1 K). In all the simulations, the cen- 
tral fragment (T > 500 K) corresponds to an adiabatic region, 
where the gas cannot radiate away its compressional energy. In 
the barotropic case, the mass of the central fragment is 4.9 x 10"^ 
Mq. The outer parts of the disk formed are much denser and 
warmer than for the two other cases, even though the gas is op- 
tically thin and should cool off" efficiently, as illustrated in the 
FLD-LF calculations. The corresponding values of the Toomre 
parameter and Jeans length are larger than in the FLD-LF case. 
Such heating, however, is spurious and reflects the approximate 
treatment of radiative cooling with a barotropic approximation, 
where temperature is set up by the density. 

The most interesting case is the FLD case with the HLLD 
Riemann solver. No fragmentation occurs in that case. The mat- 
ter falls onto a central core of mass ~ 7.3 x 10"^ Mq. With this 
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less diffusive solver, the generation of a toroidal magnetic field is 
more efficient. The interplay of magnetic field lines and velocity 
gradients leads to the effective expansion of a magnetic "tower" 
in the vertical direction (Hennebelle & Fromang 2008 ). The first 
consequence of this more robust (less diffusive) solver is that the 
disc formed is less massive in the equatorial plane than with the 
LF solver. 

Figure |2] shows the density, temperature and poloidal and 
toroidal magnetic field profiles in the z-direction at a distance 
of 10 AU from the center, for the three aforementioned calcu- 
lations, at the same time. We have noted that the extent of the 
magnetic tower depends on (increases with) the resolution. The 
magnetic tower is less dense close to the equatorial plane, since 
the gas is spread in the bubble. Contrary to the barotropic case, 
where the temperature falls to 10 K at a distance > 10 AU, the 
radiative feedback extents vertically up to ~ 100 AU with radia- 
tive transfer (FLD). Although toroidal magnetic field generation 
is more effective with HLLD than with LF, the toroidal magnetic 
field in the FLD-LF case is larger at small radii. This is due to the 
denser disk formed, which compresses and twists the magnetic 
lines in the region close to the equatorial plane. At a distance 
of 100 AU, the toroidal magnetic field component in the verti- 
cal direction is one order of magnitude larger with HLLD than 
with LF. This favors the extension of the magnetic tower. The 
poloidal magnetic field component is larger everywhere in the 
HLLD case, since matter is less compressed. 

Figure [3] displays the profiles of the magnetic braking at 10 
AU above the equatorial plane and in the z-direction at a distance 
of 10 AU from the center, for the three same calculations. In 
the equatorial plane, the magnetic braking obtained with the LF 
solver or in the barotropic case is barely significant. Magnetic 
braking favors faster accretion on the central object and thus 
yields a larger amount of kinetic energy to be radiated away at 
the first core accretion shock [all the in-falling gas kinetic en- 
ergy is radiated away at the first Larson core accretion shock). 
The gas is then significantly heated up around the central core. 
In opposite, there is more angular momentum in the (diffusive) 
FLD-LF case, which promotes fragmentation at the centrifugal 
barrier. The accretion in the central fragment is thus smaller and 
the heating due to accretion is less important. Magnetic braking 
in the z-direction is similar for the three calculations. We thus 
identify two important processes, both quantitatively affected 
by the numerical treatment: the presence of a magnetic tower, 
whose extension depends on the numerical resolution, and the 
strong magnetic braking in the equatorial plane, whose strength 
depends on the diffusivity of the solver. 

Figure |4] shows the temperature - density distribution for the 
FLD calculations with the LF solver. The isothermal and adia- 
batic regimes are recovered at low and high density (p < 10"'^ 
g cm"-' and p > 10"'^'^ g cm"-'). In between, we observe a dis- 
persion in the (T, p) plane, where matter can be hotter at low 
density and cooler at high density compared with the barotropic 
treatment (red solid line). As seen in the figure, the fragments 
eventually reach different isentropes. The central fragment lies 
on the highest entropy level, while the orbiting fragments end 
up on cooler isentropes. The barotropic EOS cannot reproduce 
such a spread in temperature at a given density. The entropy and 
the minimum Jeans mass in that case are set up by the choice 
of Pad. In the FLD-LF case, cooling is more efficient; material 
in the equatorial plane cools down by radiating in the vertical 
directions, where the gas is optically thin. However, this result is 
mainly affected by the diffusivity of the Riemann solver. In the 
FLD-HLLD case, all the material has been heated up by the cen- 
tral fragment accretion luminosity and is distributed above the 
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Fig. 2. Case yu = 20: profiles of density (a), temperature (b), 
poloidal (c) and toroidal (d) magnetic fields in the z-direction 
at 10 AU from the center, for the same calculations and at the 
same time as in fig. [T] 
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Fig. 3. Case fi = 20: profiles of magnetic breaking at 10 AU 
above the equatorial plane (left) and in the z-direction (right) at 
a distance of 10 AU from the center, for the same calculations 
(same legend) and at the same time as in fig. [T] 
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Fig. 4. Case /i = 20: temperature-density plot for each cell in the 
FLD-LF calculations (black) and with the barotropic EOS (red 
line). 



barotropic adiabat. In the FLD-LF case, the fragmentation is the 
combined result of the inefficient magnetic braking and of the 
efficient cooling due to the FLD. This illustrate the great impor- 
tance and the complexity of the interplay between magnetic field 
and radiative transfer. 



3.2. Case fi = 5 

We now consider a strongly magnetized core, with fi - 5, 
and explore the impact of radiative transfer on the outflow and 
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Fig. 5. Case // = 5: temperature (top) and density (bottom) maps 
in the xy-plane for FLD and barotropic EOS calculations (right). 



on the temperature distribution in that case. Calculations have 
been performed with the HLLD solver and Nj = 10 . For this 
more magnetized model, the magnetic field lines determine the 
dynamics of the collapse and, contrary to the previous case, it is 
the flow which ajusts to the magnetic field topology rather than 
the opposite. 

Figure 5] compares temperature and density maps in the yz- 
plane for the FLD and barotropic calculations, at f = 1.14fff ~ 
38 kyr. In both cases, an outflow is launched, with similar prop- 
agation and opening angles. The density patterns are also simi- 
lar. The temperature distributions, however, differ drastically. As 
shown in [Hennebelle & Fromang ( 2008| l, a pseudo warm disk 
is formed in the equatorial plane in the barotropic case. Note 
also the lack of significant heating in the vertical direction with 
the barotropic EOS. Since the density within the outflow is low, 
p ~ 5 X 10"'^ g cm"-', this implies a low temperature, which 
quickly falls to 11 K at 10 AU from the center. When includ- 
ing a more proper treatment of radiation, with the FLD, the ra- 
diation escapes preferentially in the vertical direction (along the 
z-axis) and heats up the gas up to about 100 AU. The highest out- 
flow velocity, ~ 1.4 km s"', is obtained in the barotropic case 
while Vr ~ 1.1 km s ' in the FLD calculation. The border line 
between the outflow and the in-falling gas shows a small spike, 
where the material is shocked and heated. The toroidal magnetic 
field component profiles are in good agreement with the results 
of Hennebelle & Fromang (2008 ) and exhibit a nearly constant 
plateau in the region of the outflow. For our present study of 
a 1 Mq core, radiative cooling in the /U = 5 case is found to 
be almost inconsequential on the collapse. The cores are never 
found to fragment. For more massive cores with stronger radia- 
tive feedback, however, the heating in the outflow is likely to be 
more dynamically significant ( Krumholz et al.|2007 1. 



4. Summary and Discussion 

In the present Letter, we have explored the effects of both ra- 
diative transfer and magnetic fields (in the limit of ideal MHD) 
on the fragmentation of a 1 prestellar core. In agreement 
with the past non-magnetic studies of, for instance, |Boss et aT] 



( |2000l i; |miFehouse & Bate] ( |2006l l; [Krumholz et al.| ( |2007] l, we 
show that a proper treatment of radiation is important to cor- 
rectly describe this process. Radiative transfer enables the gas to 
significantly cool or heat up in different regions of equal den- 
sities whereas a barotropic EOS approximation implies that the 
cooling and the heating are fixed by the density. 

In the case fj. = 20, where both the behaviour of the flow 
and of the magnetic field affect the dynamic of the collapse, we 
show that radiative transfer has an important impact on the final 
structure. A barotropic approximation cannot account either for 
the cooling of the dense and rotating gas in the equatorial plane 
nor for the heating of the less dense gas in the vertical direc- 
tion, where radiation is found to escape preferentially. Although 
a value = 20 is not favored by observations, this case clearly 
illustrates the impact of a proper treatment of radiative transfer 
on the collapse and fragmentation of prestellar cores. For the 
strongly magnetized n - 5 case, the dynamics of the collapse 
is dominated by the magnetic field. In that case, a proper treat- 
ment of radiative cooling is less consequential for the collapse it- 
self but is crucial to derive correct optical depth and temperature 
distributions, two accessible observational diagnostics. Indeed, 
the temperature distribution with the FLD is found to strongly 
depend on the geometry whereas the barotropic approximation 
yields a nearly uniform distribution. Note also that in the /i = 5 
case, we always obtain an outflow, even when using the more 
diffusive LF Riemann solver. Such outflows can not be obtained 
presently with SPH calculations. We also performed FLD test 
calculations including the radiative feedback from a pseudo cen- 
tral newborn star and found that, for our 1 Mq core case, it is 
inconsequential at the length and time scales of interest. 

The spurious diffusivity due to either the numerical resolu- 
tion or the hydrodynamic solver, however, is found to signifi- 
cantly affect the collapse and the fragmentation. A less diffusive 
numerical scheme (HLLD instead of LF) produces strong mag- 
netic braking, that transports angular momentum. Material then 
falls on the central core, leading to a larger accretion luminosity. 
A good numerical resolution is also crucial in the vertical di- 
rection to accurately describe the growth of the magnetic tower, 
which spreads the gas around the central object. 
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Appendix A: Note on the influence of tlie Riemann 
solver and of the numerical resolution 
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Fig. A.l. Case ji - 20: density maps in the xy-plane at time t ~ 
1.16 fff for 4 barotropic EOS calculations with HLLD, Nj = 6 
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In this appendix, we present a short convergence study to in- 
vestigate the effect of both, the Riemann solver and the numeri- 
cal resolution, on the fragmentation process. We show in section 
13.1! that the use of the LF Riemann solver leads to inaccurate 
results. We perform calculations using the barotropic EOS, and 
the HLLD or the LF solvers, at various numerical resolution. The 
latter is detemined by 2 parameters: A^j, the number of points per 
Jeans length, and A^exp, that gives the extent of the mesh around a 
refined cell (if a cell is flagged for refinement, then the A^exp cells 
around the flagged cell, in each direction, will also be refined). A 
large A^exp gives a smoother transition between the levels of the 
AMR grid. All calculations presented in section |3.1| have been 



run with A^exp - 4. We report results of 4 calculations, run with: 
= 1; LF, A^j = 20 and A^exp = 4; LF 



HLLD, A^j = 6 and A^e: 
A^j = 10 and A^exp = 1; LF, A^j = 6 and A^exp = L 



Figure A. 1 portrays the density maps in the equatorial plane 
for the 4 calculations, at time t ~ 1.16ffr. The calculations per- 
formed with the LF solver, Nj - 6 and A^exp - 1 clearly diverge 
from the other ones, since it fragments. The lack of resolution 
clearly induces inaccurate fragmentation. The three other calcu- 
lations are qualitatively similar (no fragmentation), and in accor- 
dance with Hennebelle & Teyssier ( 2008[ ) results (obtained using 
a Roe type Riemann solver) . This indicates that the numerical 



resolution has to be enhanced with the LF solver in order to avoid 
spurious eff'ects due to the diffusivity of the solver. 

Eventually, note that these differences are more important in 
the case of FLD calculations, since radiative transfer can have 
two opposite effects, depending on the strength of the magnetic 
braking. If the magnetic braking is not significant (with a bad 
resolution and a diffuse Riemann solver), the material at the 
centrifugal barrier cools down and then fragments. On the other 
hand, when the magnetic breaking is efficient (with HLLD or a 
high resolution), the infall velocity and consequently, the accre- 
tion luminosity, become larger, which prevents fragmentation to 
occur (material looses angular momentum and heats up). 

Appendix B: Case /u - 1000 

For this quasi-hydro case, the initial dense core is highly 
gravitationally unstable. We have performed three types of cal- 
culations with the LF Riemann solver: one with the barotropic 
EOS and y = 5/3, and two with the FLD and diff'erent adia- 
batic exponents (y - 5/3 and y - 7/5). The value y = 7/5 is 
more appropriate when T > 100 K, when the rotational degrees 
of freedom of H2 are excited (e.g. Machida et al.|2008 i. Using 
y - 7/5 yields cooler adiabatic cores than with y = 5/3 at the 
same density. Higher densities and temperatures in the core are 
thus reached more rapidly. 



Figure A. 2 portrays density and temperature maps in the 
equatorial plane for the three aforementioned calculations at a 
time f - 1.15 fff. Temperature maps range from 1 1.4 K to 20 K. 
As pointed out in Commergon et al. ( 2008 ), the horizon of pre- 
dictability is very short when the mitial thermal support is low. 
Calculations do not show a convergence to the same fragmen- 
tation pattern, since we do not integrate the same system equa- 
tions. The FLD calculations yield different fragmentation modes 
for the two different y's: with y - 7/5, we get a central object 
and two satellites (4 at later times), whereas with y = 5/3, we get 
a ring fragmentation pattern, with orbiting fragments linked by a 
bar. The barotropic EOS calculation does not produce a central 
object. From the temperature map of the FLD case with y = 7/5, 
we see that the central region is quite hot. Fragments are more 
compressed and they radiate more energy. Fragments are also 
formed faster than the time required for the radiative feedback 
to become significant (beginning of the second collapse). 
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Fig.A.2. Density and temperature maps in the equatorial plane at time t= 38.3 kyr (l.lSfff), for the 3 unmagnetized cases: FLD 
approximation and 7 = 7/5 (left), the FLD approximation and y - 5/3 (middle) and with a barotropic EOS using p^d - 2.3 x 10"'^ 
g cm"^ and y = 5/3 (right). Velocity vectors are over-plotted on the density maps 



